Short-ranged attractions in jammed liquids: How cooling can melt a glass 
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We demonstrate that an extended picture of kinetic constraints in glass-forming liquids is suffi- 
cient to explain dynamic anomalies observed in dense suspensions of strongly attracting colloidal 
particles. We augment a simple model of heterogeneous relaxation with static attractions between 
facilitating excitations, in a way that mimics the structural effect of short-ranged interparticle at- 
tractions. The resulting spatial correlations among facilitated and unfacilitated regions give rise 
to new relaxation mechanisms that account for non-monotonic dependence of relaxation times on 
attraction strength as well as logarithmic decay of density correlations in time. These unusual fea- 
tures are a simple consequence of spatial segregation of kinetic constraints, suggesting an alternative 
physical perspective on attractive colloids than that suggested by mode-coupling theory. Based on 
the behavior of our model, we predict a crossover from super-Arrhenius to Arrhenius temperature 
dependence as attractions become dominant at fixed packing fraction. 



Homogeneous liquids with slowly varying intermolec- 
ular attractions are typically well represented by ap- 
propriate hard sphere reference systems [lj. The deeply 
supercooled regime of simple liquids appears to be no 
exceptional. Computer simulations of hard sphere fluids 
at high packing fraction, cj>, exhibit the canonical fea- 
tures of glass-forming materials, including dynamic het- 
erogeneity, stretched exponential decay of dynamical cor- 
relators, and dramatic changes in relaxation times fol- 
lowing small changes in density (correspondingin most 
experiments to small changes in temperature T) 3] . Con- 
focal microscopy of colloids under conditions designed to 
maximally screen interparticle attractions confirms this 
fact0. Most of the theoretical approaches developed to 
rationalize sluggish dynamics take advantage of the fact 
that intermolecular structure at such high densities (or 
low temperatures) differs little from that of the standard 
liquid state near the triple point. It is exactly this spa- 
tially uniform structure that justifies a hard sphere rep- 
resentation. 

Glass-forming liquids that do not fit into this van der 
Waals picture thus not only introduce the possibility 
of qualitatively new dynamical features but also pro- 
vide a significant test of the flexibility of various the- 
oretical perspectives. The simplest such material is a 
fluid of hard spheres that attract one another strongly 
over short distances, realized experimentally as a sus- 
pension of attractive colloidal particles. In this case 
attractive forces can play an important role in deter- 
mining intermolecular structure. Experiments and com- 
puter simulations indicate that for small u/T (where u is 
the strength of attractions), overall relaxation rates ex- 
ceed those of a hard sphere liquid at the same packing 
fraction!! 0, 0, El EH EH • For fixed u, one can thus 
view melting of the glass as a consequence of cooling. For 



large u/T, however, tight clustering of particles leads to 
rigid gel-like structures that rearrange extremely slowly. 
This second, attraction-driven glassy state would seem 
to differ qualitatively from the original jammed material 
(u = 0) in both structure and dynamics [12J. 

Calculations based on mode-coupling theory have mo- 
tivated much of the experimental study of slow dynam- 
ics in attractive colloids, predicting the existence of mul- 
tiple glassy states and relaxation times that vary non- 
monotonically with u/T at fixed ^0 El El- The 
physical origins of these predictions, however, are dif- 
ficult to assess. Attraction-driven changes in dynam- 
ics have been ascribed to loosening of the dense liq- 
uid environment that confines particle motion. But the 
mean field approximation underlying mode-coupling the- 
ory would seem unable to capture the inhomogeneous 
particle clustering necessary for such a mechanism. This 
mean field nature is most problematic in generating a 
spurious dynamical arrest at large u/T or <f>. Nonethe- 
less, many predictions of mode-coupling theory, includ- 
ing logarithmically slow relaxation near re-entrance to 
the attractive glass, have been verified by dynamic light 
scattering H 00, H ©, El EH 

A more intuitive, though less detailed, theory of deeply 
supercooled liquids has been constructed by considering 
only the inhomogeneous fluctuations that lie outside the 
scope of mode-coupling theory [lH EH Eil . Jamming 
at high density prevents all but a small fraction of par- 
ticles from moving a significant distance over any short 
time interval. The loose regions (or defects) which enable 
nearby relaxation are therefore sparse and essentially un- 
corrected in space. Although the statistics and spatial 
patterns of these defects may be trivial, their dynamics 
are complicated by the fact that jammed regions can- 
not spontaneously loosen unless nearby motion provides 



2 



space for them to do so. In other words, defects can 
only be created in the vicinity of other defects. Such a 
constraint greatly limits the accessible pathways for ex- 
ploring phase space . 

The simplest model embodying this picture consists of 
non-interacting excitations whose dynamics are coupled 
together by kinetic constraints. Each cell of a lattice, 
with the size of a liquid's density correlation length, is 
thus designated as either excited, rii — 1, or jammed, 
rii = 0. Creating defects in a liquid requires local reor- 
ganization of particles and thus costs entropy (or, in a 
thermal system, free energy). This cost is represented 
schematically by an external field, h > 0, limiting the 
average concentration of mobile cells, (rii) = c <C 1. 
The model system's total energy is then TL — h^iLi n i- 
(Throughout this paper we express all energies in units of 
T, distances in units of the lattice spacing, and frequen- 
cies in units of the rate at which an unbiased, uncon- 
strained cell changes its excitation state.) Cell mobilities 
change with rates that preserve a canonical distribution, 
but are subject to a constraint of facilitation. Namely, 
cell i is only free to change its state when one of its near- 
est neighbors is already excited. In the one-dimensional 
East model only cell i — 1 can facilitate changes at cell 
i |2f]j| . In the Fredrickson- Andersen model either neighbor 
(i + 1 or i — 1) can provide facilitation[2l|. 

Despite their simplicity one-dimensional facilitation 
models exhibit remarkably nontrivial dynamics. Garra- 
han and Chandler have shown that kinetic constraints 
in these models confer special geometries on trajecto- 
ries below a crossover value of c0|. These spatio- 
temporal patterns result in sluggish relaxation with tem- 
perature dependences mirroring those of many glass- 
forming materials^JH^,!^. There are even indications 
that scaling of such geometric features near a dynamical 
critical point at T — may be universal among facilita- 
tion models and real liauids[2^|. Nevertheless, a detailed 
relationship between between mobility-promoting defects 
and molecular degrees of freedom remains elusive. 

It is not obvious that the kinetic facilitation picture 
is sufficiently flexible to account for the anomalous dy- 
namical behavior reported for colloidal suspensions |25|. 
Basic facilitation models possess a single control param- 
eter, c, which seems best to correspond to the density (or 
temperature) of a real material. The monotonic growth 
of relaxation times with decreasing c in these models 
appears to prohibit a straightforward explanation of re- 
entrance and multiple glassy states in terms of dynamic 
heterogeneities alone. In this Letter we describe a sim- 
ple extension of the kinetic facilitation picture inspired 
by structural changes arising from short-ranged attrac- 
tions. We then examine the dynamics of an appropriately 
modified one-dimensional East model and demonstrate 
striking agreement with the phenomenology of attractive 
colloids. 

We assert that the primary effect of introducing short- 



ranged attractions in a dense, disordered material is to 
produce a gradual segregation of tightly clustered, im- 
mobile regions from loose regions enriched with facili- 
tating defects. In the limit of very strong attractions 
such demixing produces gel-like structures and eventu- 
ally macroscopic phase separation, as seems appropriate 
for strongly attractive colloids. In the context of one- 
dimensional facilitation models, gradual demixing is nat- 
urally generated by adding attractions between adjacent 
excitations: 



N 



H = -e 2J nin. i+ i + h(e, c) 2J rii 



(1) 



Although clusters and voids would certainly segregate as 
desired with increasing attraction strength e, the aver- 
age concentration of mobile cells would increase signifi- 
cantly if the external field h did not change. We imag- 
ine that increasing attraction strength in a liquid with 
fixed density, however, does not significantly change the 
concentration of facilitating defects (in a sense, the to- 
tal free volume). If we view c as uniquely determined 
by 4>, then the field strength h must be adjusted to 
maintain a desired value of c for any e. With a routine 
change of variables, Sj = 2rii — 1, Eq. ^ becomes (within 
an additive constant) the energy of a one-dimensional 
Ising model with exchange interaction J = e/4, mag- 
netic field H = —(h — e)/2, and magnetization per spin 
(si) = — l + 2c. The required adjustment of external field 
given c and e is thus a standard result |26| : 



h(e, c) = e + 2 sinh 



1 - 2c 



e/2 



(2) 



This constraint on the density of excitations, motivated 
by the notion of conserved free volume at a specific 
packing fraction, distinguishes our approach from sim- 
ply adding attractions to a conventional model of kinetic 
facilitation. It is an essential feature, without which dy- 
namical anomalies such as re-entrance would not be ob- 
tained. 

We emphasize that the attraction between defects in 
our model does not represent a literal interaction energy 
between mobile regions in a liquid. In a liquid it is of 
course particles that attract one another, the effect of 
which is to drive mobile regions together. Attractions in 
Eq. could thus be viewed as a potential of mean force 
between mobile regions, or simply as a measure of cluster- 
ing tendency. For a fixed value of c, it is in fact energet- 
ically immaterial whether attractions act between adja- 
cent excited lattice cells, rii, or between adjacent jammed 
cells, I — rii. 

The structural effect of nonzero e for a particular value 
of c is clear: Excitations cluster together, and the average 
spacing between consecutive clusters, d, grows: 



d + (e £ - 1) 

c -l/2 e e/2 



e < h 
e » h . 



(3) 
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FIG. 1: A representative trajectory of the one-dimensional 
attractive East model with c = 0.10 and e = 4.0. Vertical 
cross-sections are snapshots of a portion of the system (sev- 
eral hundred lattice cells), in which excited lattice cells are 
colored black. Unexcited cells are not shown. Time runs 
along the horizontal axis. The total duration of this trajec- 
tory segment is 1000 units of the relaxation time of a single 
facilitated cell. See Ref. [lSJ for a detailed discussion of this 
manner of depicting trajectories and the space-time geometry 
it reveals. 



Here, do ~ c _1 and ho « ln(l/c) are the values of d and h 
at e = 0. Just as short-ranged attractions enhance tran- 
sient structural heterogeneity in a liquid, attractions in 
the East model generate nontrivial (i.e., non- ideal gas) 
spatial patterns of defects. The dynamical consequences 
are more subtle. In a Metropolis Monte Carlo trajectory 
the rate of exciting a cell at the boundary of an excita- 
tion domain (or alongside a single excitation) is e~( h ~ e ^ . 
Larger values of J in the isomorphic Ising model require 
weaker symmetry-breaking magnetic fields to maintain 
fixed magnetization, so that h(e, c) — e decreases mono- 
tonically with e: 



h - e + 2c(e c - 1), e < h (A , 



Attractions therefore expedite fluctuations in the vicinity 
of excitations. This enhanced defect creation rate must 
be balanced, however, in order to preserve the proper 
equilibrium distribution. Here, balance is provided by 
the decline in number of excited clusters per lattice cell. 

In the FA model defects may diffuse freely through 
the system through short sequences of local moves, e.g., 
10 — ► 11 — ► 01. The directionality of constraints in the 
East model gives rise to more complex, self-similar re- 
laxation pathways when the density of excitations is low. 
In this case all excitation domains are pinned on one 
side, as is clear from the trajectory plotted in Fig. 1, 
and can fully relax only when contacted by another do- 
main on that side. For e = 0, the minimum energy 
path for connecting two defects separated by a distance 
I = 2™ establishes a set of isolated defects between 
them. These stepping stones are situated hierarchically, 
at f/2, 31/4, 71/8, . . . , £ — 1, giving rise to an energetic 
barrier A(^) = mho 27]. The overall relaxation rate is 
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FIG. 2: Renormalization of the hierarchical relaxation mech- 
anism characteristic of the East model. A small portion of 
the system is depicted horizontally at each step. Numbers (0 
and 1) indicate the presence or absence of individual excita- 
tions in a sequence of configurations proceeding from top to 
bottom. Each of the 10 steps comprising this schematic tra- 
jectory involves creating or (nearly) destroying blocks of a = 8 
consecutive excitations. (The simple sequences of single-cell 
changes that create and destroy these blocks are not shown.) 
This idealized example represents the lowest energy pathway 
that results in elimination of the rightmost defect. 



e -A(do) = di 



-h /ln2 



then To ~ c • • — u 

Attractions establish a length scale, a « hj(h — e), be- 
low which it is energetically advantageous to fill immo- 
bile spaces between excitations with additional defects. 
Gap filling (clearly visible in the trajectory of Fig. 1) 
effectively renormalizes the basic hierarchical relaxation 
mechanism. Low energy pathways, such as that depicted 
in Fig. 2 for a — 8, do not involve isolated defects sepa- 
rated by distances smaller than a. This change of scale, 
which grows monotonically with e, 



l + e/ho, 
ecVV/ 2 , 



e < h 
e » h , 



(5) 



has two important consequences. First, the activation 
energy for pathways connecting subsequent excitations 
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renormalizes. For £ — (2 m — l)a + 1, the barrier be- 
comes A(£) = (m-l)h + a(h-e) ps h\n{£/a)/\n2. Sec- 
ond, a new time scale r a emerges, corresponding to the 
frequency with which excitation domains spontaneously 
grow to size a. 

Because in the East model domains are bound on one 
side, fluctuations of the boundaries dividing excited do- 
mains from immobile domains (i.e., domain walls) are 
equivalent to the motion of a one-dimensional random 
walker. The potential energy of this walker at position x 
is (h — e)x, and a reflecting boundary restricts motion to 
x > 0. Excursions for which the change in energy is com- 
parable to T are essentially diffusive in nature, occurring 
with a frequency proportional to (h — e) 2 . Excursions 
of length L > (h — e)^ 1 are less frequent by a factor of 
exp[(/i — e)L — 1]. The basic time scale associated with 
steps of length a ;§> 1 is therefore r a ps (h — e)~ 2 e h . 

Assembling these results, we estimate the relaxation 
rate for e > as 



r w o 



\2 

e) exp 



In 2 



ln(d/a) 



(6) 



This expression depends on attraction strength e implic- 
itly through h, d, and a. Each of these quantities grows 
with increasing e, but their combination in Eq.[^]may be 
nonmonotonic. 

Through rcnormalization of the basic dynamical length 
scale, attractions aid the propagation of excitations. But 
for fixed c, attractions also increase the distance excita- 
tions must propagate to achieve relaxation. The quantity 
h\n{d/a) characterizes this competition. For weak at- 
tractions, ln(d/a) ps Indo—e/tiQ varies much more rapidly 
than h ~ ho + 2ce, and relaxation becomes more facile 
with increasing e: 



In [ — 

To 



mT £<<V 



(7) 



Defect concentration c does not appear in this result, 
which should be accurate in the limit of small c. Terms 
of order (lnc) -1 , which for moderately small values of 
c cause the initial growth in lnT to differ for different 
attraction strengths, have been omitted. 

For strong attractions, ln(d/a) ps ln(e/c) varies instead 
more slowly than h ps e, so that T decreases with increas- 
ing e. In this regime the change of dynamic length scale is 
insufficient to offset the enhanced sparseness of excitation 
domains. Relaxation rates should thus be maximum at 
the crossover from small-e to large-e behavior, i.e., near 
e = ho- 

The hierarchical mechanism we have described is only 
sensible for a < d. When the renormalized length scale 
is comparable to or larger than d, there are few con- 
secutive defects separated by distances larger than a. 
For e ^ c~ x the minimum energy pathway connecting de- 
fects simply excites all intervening lattice cells. Domain 




FIG. 3: Average time r required for all but N/e 
cells to flip, as a function of attraction strength e. 
The concentrations are, from bottom to top, c = 
0.02,0.04,0.06,0.08,0.1,0.12,0.14,0.16, and 0.18. The con- 
tinuous time Monte Carlo simulations that generated these 
results employed TV = 10 lattice cells. 



wall fluctuations are then the only relevant dynamics. 
A domain wall excursion of length d requires activation 
energy d(h — e), which approaches c _1 in the limit of 
large e. Relaxation rates for a ^ d vary with attrac- 
tion strength (at fixed c) in a simple Arrhenius fashion, 
r« [{h-e)] 2 e 1 / c ^c~ 1 e^ +1 ' c . 

Numerical simulations of the East model defined by 
Eqs.n an d 12] verify these expected features. Relaxation 
time r is plotted in Fig. 3 as a function of e for several 
values of c between 0.01 and 0.2. We define r through 
the persistence function 7r(t), i.e., the average fraction of 
cells that do not change over a time interval of length 
t[23|. Specifically, 7t(t) = 1/e, so a majority of cells have 
undergone fluctuations within time r. Most importantly, 
numerical results confirm the existence of re-entrant dy- 
namics as attraction strength increases. The maximum 
relaxation rate can be nearly two orders of magnitude 
larger than its value at e = in this range of defect 
concentration. As expected, this maximum occurs near 
e = ln(l/c). Re-entrance should become still more dra- 
matic for smaller values of c. 

For the larger values of c we have considered, the 
asymptotic behavior for a ^ d begins to appear in our 
simulation results at large e. (For smaller c, the corre- 
sponding range of e is numerically inaccessible.) Several 
results point to dominance of activated domain wall fluc- 
tuations in this regime. First, the Arrhenius behavior 
of lnT evident in Fig. 3 is consistent with our analy- 
sis. Interestingly, although the hard sphere liquid (and 
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the e = East model) is a fragile glass-former |28j . our 
results suggest that the re-entrant glass state is instead 
strong. This prediction seems reasonable, given that gel- 
like structures have much in common with convention- 
ally strong, network-forming materials [2!|. It would be 
straightforward to test in experiments or molecular dy- 
namics simulations by varying temperature at fixed pack- 
ing fraction and attraction strength. 

The decay of the persistence function exhibits a more 
direct signature of domain wall fluctuations. For small 
values of c, a great majority of lattice cells are unex- 
cited at any time. The average fraction of relaxed spins, 
1— 7r(i), should thus be largely determined by the growth 
of excited domains. (Relaxation of the initially excited 
minority contributes only weakly and at early times.) Re- 
call that the time required for a typical domain wall ex- 
cursion of length L is t « (h — e)~ 2 e( h ~^ L . In other 
words the typical extent of domain growth at time t is 
L ~ (h — e) _1 \n[(h — e) 2 t]. As a result, the persistence 
function, 



1 - L/d 
constant 



(8) 

cln(t) (9) 

decays logarithmically. Eq. [5] suggests that the coeffi- 
cient of logarithmic relaxation for n(t) and related cor- 
relation functions should not change significantly with 
u/T for fixed packing fraction. ir(t) is plotted in Fig. 4 
for c ~ 0.14 and several values of e. A period of loga- 
rithmic decay, with a slope on the order of c, appears as 
attractions are introduced and becomes more prominent 
as e increases. Since the average spacing between excited 
domains still exceeds a in this case, the terminal decay 
of 7T (i) maintains a form characteristic of hierarchical re- 
laxation. Logarithmic decay should extend over longer 
time spans for smaller c and larger e. The two decades 
of such decay apparent in Fig. 4 for e > 6 reflect only the 
limit of our computational resources. 

In summary, the simplest physically motivated modi- 
fication of the East model reproduces the principal dy- 
namical anomalies of attractive colloidal suspensions. In 
contrast to mode coupling theory, our model suggests 
simple physical origins of these features. It is the com- 
bination of static and dynamic heterogeneity in our pic- 
ture that generates new dynamical features. Clustering 
of mobility renormalizes the relaxation mechanism that 
predominates in the absence of attraction. Re-entrance 
results directly from a competition between a growing 
dynamic length scale a and a growing static length scale 
d. When the scale of uninterrupted defect propagation 
exceeds the typical separation between defect domains, 
the primary mode of relaxation crosses over to simple do- 
main wall fluctuations, whose dynamics give rise to log- 
arithmic relaxation. Although we have presented results 
only for a modified East model, an equivalcntly modified 
FA model shares much of the same behavior. Basic re- 
laxation pathways are similarly renormalized, leading to 
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FIG. 4: Persistence n(t) as a function of time t for several 
values of e with c = 0.14. Dashed lines highlight regions 
where logarithmic relaxation is noticeable. 



re-entrance of comparable magnitude. Logarithmic re- 
laxation is less pronounced in this case, since excitation 
domains are not pinned on one side. At long times, fluc- 
tuations of a domain's two walls become correlated, and 
our arguments for the decay of 7r(i) do not strictly apply. 

Based on the dynamics of the attractive East model 
we make several predictions that may be tested by ex- 
periment or by numerical simulation of colloid dynamics. 
Most significantly, our results imply for increasing attrac- 
tion strength not only a transformation from strong to 
fragile glass but also rapid growth of the length scale as- 
sociated with dynamic heterogeneity. Recent computer 
simulations 30, 131| are consistent with this latter predic- 
tion. 

We are aware of another effort to rationalize the be- 
havior of attractive colloids within the kinetic facilita- 
tion picture 32] . This approach, in contrast to our own, 
preserves a lack of static spatial correlations between de- 
fects. It achieves re-entrance by instead modulating the 
details of facilitation constraints. The two pictures are 
most distinct for very strong attractions, where spatial 
correlations between defects in our model become long- 
ranged. In our view this behavior naturally reflects the 
growing range of correlated density fluctuations in a dis- 
ordered material as the strength of interparticle attrac- 
tions increases (or temperature decreases). Detailed ex- 
periments or computer simulations should be sufficient to 
assess the relative validity of these approaches by charac- 
terizing spatial distributions of dynamic heterogeneities 
on short time scales. 

Although the kinetic facilitation picture is most sen- 
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sible at high density, it is interesting to consider the 
implications of our results for low packing fraction and 
very large attraction strength. In a version of the at- 
tractive East model in more than one dimension, nearly 
irreversible aggregation of immobile regions would de- 
scribe the growth of sparse fractal patterns similar to 
those observed in experiments |33l l34|. It may therefore 
be possible to map the entire phase diagram of our model 
(in the plane of c _1 and e) onto that of a real attractive 
colloidal suspension (in the plane of </> and u/T). 
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